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To capture various experimental results in the pseudogap regime of the underdoped cuprate 
superconductors for temperature T < T*, we propose a four-component pair density wave (PDW) 
state, in which all components compete with each other. Without random field disorders (RFD), 
only one of the PDW components survives. If the RFD is included, this state could become phase 
separated and consist of short range PDW stripes, in which two PDW components coexist but differ 
in magnitudes, resulting in charge density waves (CDW) and a time-reversal symmetry breaking 
order, in the form of loop current, as secondary composite orders. We call this phase-separated pair 
nematic (PSPN) state, which could be responsible for the pseudogap. Using a phenomenological 
Ginzburg-Landau approach and Monte Carlo simulations, we found that in this state, RFD induces 
short range static CDW with phase-separated patterns in the directional components and the static 
CDW is destroyed by thermal phase fluctuations at a crossover temperature Tco < T* , above which 
the CDW becomes dynamically fluctuating. The experimentally found CDW with predominantly 
d-wave form factor constrains the PDW components to have s' ± id pairing symmetries. We also 
construct a lattice model and compute the spectral functions for the PSPN state and find good 
agreement with ARPES results. 

PACS numbers: 71.45.Lr, 74.20.Rp, 74.72.Kf, 74.81.-g 


I. INTRODUCTION 

Since the discovery of the cuprate superconductors [1] 
(SC), the pseudogap regime has attracted a lot of atten¬ 
tions for its exotic properties [2] and the possible rela¬ 
tion to the high temperature d-wave superconductivity 
(dSC). Two main types of theories are proposed for the 
pseudogap. One type [3-6] suggests that pseudogap is a 
precursor to the dSC phase with pre-formed Cooper pairs 
without global phase coherence. However, mounting ex¬ 
perimental evidences are pointing otherwise that pseudo¬ 
gap is indeed another broken symmetry state. The recent 
inputs mainly come from the X-ray scattering experi¬ 
ments [7-15], which suggest the state is characterized by 
the onset of incommensurate charge density wave (CDW) 
with wave vectors at (±2(5, 0) and (0, ±2(5), where Q de¬ 
creases with increasing hole doping. This finding is in ac¬ 
cord with earlier scanning tunneling microscope [16-20] 
(STM) and nuclear magnetic resonance [21, 22] (NMR) 
experiments. There are some early theoretical studies to 
explain the pseudogap in terms of CDW [23, 24]. The 
phase diagram is however much richer. Other orders as¬ 
sociated with different broken symmetries are also de¬ 
tected in the pseudogap, e.g., the time reversal symmetry 
breaking (TRSB) order detected in polar Kerr rotation 
[25, 26] (PKR) and polarized neutron diffraction [27-30] 
(PND) experiments and the nematic order [13, 31-34] 
that breaks the C 4 lattice rotational symmetry. This 
complicated zoo of orders poses a natural question that 
whether we can unify these into a common origin and un¬ 
derstand the pseudogap in a coherent manner. Among 
the various theories proposed, pair density wave (PDW) 
order, which is a spatially modulating SC state similar to 
Fulde-Ferrell-Larkin-Ovchinnikov (FFLO) states [35, 36], 
previously studied in different contexts for the cuprates 


[37-46], is suggested to be responsible for the pseudogap. 
Recent works [47, 48] by Lee and Agterberg et al. have 
successfully explained many features of the pseudogap. 
In the PDW theory, the CDW is induced by the PDW 
as a secondary composite order, accounting for the STM 
and X-ray scattering results. Moreover, the PDW can 
also induce a loop current (LC) order [48-50] to account 
for the TRSB order observed in PKR [25, 26] and PND 
[27-30] experiments. Nonetheless, the PDW order also 
explain several ARPES features [12, 51-55], namely the 
kp-kc misalignment, antinodal gap closing from below, 
and the Fermi arcs. 

Hinted by these successes, the PDW order might hold 
the key to understand the pseudogap of hole doped 
cuprates. However, several issues have surfaced. Most 
noticeably, the charge density waves observed in STM 
[19, 20, 56] and REXS [13] are short ranged. STM fur¬ 
ther show that these density waves form a domain struc¬ 
ture and are directional within each domain. The same 
conclusion is also drawn from REXS. While these results 
are obtained in no magnetic field, the high field exper¬ 
iments, on the other hand, reveal that the CDW with 
the same in-plane wave vector is long ranged and unidi¬ 
rectional [21, 57-59]. As a first step to understand this 
onset of long range CDW state at high field [22, 58, 59], 
it is essential to explain why the CDW observed in weak 
field is short ranged. 

Experiments however impose stronger constraints, be¬ 
sides the possible order parameters, on theories. Each 
set of experiments reveals a specific doping dependent 
temperature scale but not all these temperature scales 
can fall into a single simple curve in the phase diagram. 
This challenges the belief that the pseudogap can be as¬ 
cribed to a single origin. We shall mention those that 
are relevant in the current paper. The electronic trans- 
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FIG. 1. (Color online) Bare Fermi surface. The PDW pairing 
centers ±Qx,2//2 are at the intersections with the Brillouin 
zone boundary. The blue lines are for the quasiparticle spec¬ 
tral function scans in Fig. 9 and 10. 


port giving rise to linear resistivity [60, 61] defines Tp 
that is believed to capture a quantum critical point as¬ 
sociated with corresponding order parameter and it is 
commonly taken as the definition of the pseudogap tem¬ 
perature scale T* = Tp. ARPES measurements define 
^ARPES ^ below which unusual quasiparticle spec¬ 
trum is revealed. Moreover, REXS [7, 9, 11-15] also 
found the fluctuating CDW correlation below T*, but 
RIXS [8, 10] and STM [17] reveal a lower temperature 
scale Too < 7^* for static charge order. It is found that 
Tco has a maximum at some doping but T* decreases as 
doping increases. A valid theory would be able to explain 
the emergence of these temperature scales within the 
pseudogap. Similar to Tco^ we note that for TRSB or¬ 
der, there are also two distinct temperature scales found 
in PND [27-29] (see also Ref. [30]) and PKR [25, 55] (see 
also Ref. [26]) experiments. The PND experiments found 
that below Tm ^ 7"* an intra-unit cell (lUC) TRSB mag¬ 
netic order exists, while the PKR experiments, which also 
detect TRSB, reveal a strictly lower temperature scale 
Tk < Tm but doping dependence of Tk shows a similar 
trend as Tm- However, we shall address the issue of Tk 
elsewhere. 

Motivated by these issues, we study the effect of ran¬ 
dom field disorders [62, 63] (RED) in a four-component 
PDW (PDW4) model [44, 46, 47] (see also Ref. [64] for the 
study of a CDW model). Intuitively, RED does nothing 
more than inducing the PDW short range [65]. But we 
will show below that in order to explain the experimen¬ 
tal results coherently, inclusion of RED would lead us to 
consider a PDW state of different nature unexplored be¬ 
fore, thus it is not straightforward to generalize existing 
PDW results to our case. We consider four PDW compo¬ 
nents /S.±Q^ y at wave vectors and TQy and all the 
components compete with each other. The pairing cen¬ 
ters TQxl2 and TQyj2 are located at the intersections 
of the Eermi surface and the Brillouin zone boundary 
as shown in Eig. 1. We imagine that the PDW orders 



FIG. 2. Schematic phase diagram for hole doped cuprates, 
where <5 is the doping and T is temperature. T* is the char¬ 
acteristic temperature for the pseudogap, in which we pro¬ 
pose the phase-separated pair nematic state is responsible for. 
Tw is the critical temperature for the antiferromagnetic order 
near zero doping. Focusing at a doping at underdoped 
regime, Tc is the critical temperature for the dSG, Tco is for 
the presence of short range static GDW. 


present below a characteristic temperature T* that de¬ 
creases as doping increases, as shown in Eig. 2. T* termi¬ 
nates inside the dSC phase at a quantum critical point 
6 c, which controls a quantum critical region [66-68] with 
special properties like linear electronic resistivity [60] (see 
also Ref. [61]). We first use Ginzburg-Landau (GL) for¬ 
mulation and Monte Carlo simulations to study the in¬ 
terplay of PDW, CDW and RED in real space. A PDW 
with ±Q components would induce a secondary compos¬ 
ite CDW order with wave vector 2 Q [41, 44-46], which is 
then coupled to the RED [65]. With the TQ PDW com¬ 
ponents compete with each other and the presence of 
RED, it is possible to induce a TRSB loop current order 
compatible with the CDW. We also show that there ex¬ 
ists a crossover temperature scale Tco (see Eig. 2) for the 
short range static CDW order and argue that the PDW 
superconductivity is absent due to SC phase fluctuations. 
Moreover, RED induces a domain pattern, or phase sepa¬ 
ration, on the directional PDW components at TQx and 
TQy. We shall dub this state of short range PDW with 
four competing components the phase-separated pair ne¬ 
matic (PSPN) state. Next, we argue that the predomi¬ 
nantly d-wave form factor CDW observed in STM con¬ 
strains the PDW to be a bond order with s' ± id-wave 
pairing symmetry. Using this input together with a lat¬ 
tice model and the GL functional, we compute the quasi¬ 
particle spectral functions under the influence of thermal 
fluctuations and find good agreement with ARPES. 

In this study, we make several assumptions. We con¬ 
sider only the cuprate families that show CDW of de¬ 
creasing wave vectors with increasing doping, as seen in 
recent X-ray scattering experiments [10, 12, 14]. In this 
study, we have neglected the cuprate families like LSCO 
and LB CO that show CDW wave vectors scaling with 
doping [38] and the CDW form factor appears to be s'- 
wave [69], though PDW has also been suggested to play 
a role [40, 43, 44, 65]. We shall focus on the experimen¬ 
tal results in weak magnetic fields, because some recent 
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experiments [21, 22, 57-59] reveal that the CDW order 
observed in high field and low temperature could be of 
different nature than the one obtained in weak field. In 
the following, we neglect the dSC for simplicity and only 
consider the PDW order and its possible induced com¬ 
posite orders at a specific doping in the phase diagram 
(Fig. 2). It is suggested that the d-wave superconductiv¬ 
ity competes with the CDW [70-73]. Lastly, the origin 
of the PDW in cuprate is not the focus of the present 
paper, but some previous studies show that PDW could 
be originated from strong correlations [47, 74]. 

There is another line of theory that consists PDW as an 
important player. Its starting point is a spin-fluctuation 
scenario, captured within the spin-fermion model [75]. 
PDW and CDW with d-wave form factors can emerge 
naturally from this semi-microscopic model [76-79] (see 
also Ref. [80]) and it can also capture a number of impor¬ 
tant properties and features observed in the experiments. 
In particular, an unidirectional CDW/PDW mixed state 
[77, 78] obtained from the model can account for the 
ARPES spectrum, TRSB order, CDW order, and break¬ 
ing of C 4 symmetry. Although these results are similar 
to those of the PSPN state we discuss here, we note that 
in the spin-fermion model, PDW and CDW orders form 
a “supervector” stemming from an approximate SU(2) 
symmetry near the intersection points of the Fermi sur¬ 
face and magnetic Brillouin zone (the “hot-spots”), while 
here the PSPN state has the PDW as the main character 
and the CDW plays only a parasitic role [47, 48]. As we 
will show below, this parasitic CDW plays a crucial part 
in explaining the distinct temperature scales observed in 
the pseudogap, namely T* and Too- Here we wish to 
consider the parasitic CDW induced by the PDW only 
and see how far this perspective can proceed. 

The paper is outlined as follow: In section II, we will 
first briefly review the two-component PDW state, and 
followed by the discussion of the interplay between the 
four-component PDW state and RED in section III. In 
section IV, we constrain the PDW pairing symmetries 
from the CDW experiments. In section V, we present the 
results for the quasiparticle spectral functions and com¬ 
pare with the ARPES experiments. Due to the number 
of acronyms, we provide an acronym list for the phrases 
commonly used at the end of the paper to facilitate the 
reading. 


II. BRIEF REVIEW ON TWO-COMPONENT 
PDW (PDW2) 

In this section, we briefly review and discuss the model 
of two-component PDW (PDW2) and the effect of RFD. 
The main focus of the paper is to study the PDW with 
four components and relate to the experimental results 
in the pseudogap regime of the cuprates. In order to do 
that, it is beneficial to first discuss the PDW2 system 
[38, 40, 43-46, 65], from which we can easily generalize 
the results to the four component case. Firstly, we will 


introduce the PDW order parameters consisting two com¬ 
ponents and the corresponding GL functional. Also it is 
possible to induce a CDW from the PDW2, provided that 
both of the PDW components are non-vanishing. Next, 
we imagine that the CDW are coupled to the charged 
impurities, which are modeled by the coupling between 
CDW and RFD in the GL functional. For a moment, we 
shall consider only the subproblem of the CDW coupled 
to RFD, which can be cast into the random field XY 
model [81, 82]. (Actually the random field XY model 
and variants have been taken as minimal models of a 
broad class of condensed matter systems with quenched 
disorders, such as flux-line arrays in dirty type-II super¬ 
conductors [83], charge density waves [84], and smectic 
liquid crystals in random environments [85-87]) In the 
model, the XY phases (CDW smectic phase) are coupled 
to the RFD. It is found [81, 82] that in 2D the RFD, no 
matter how weak, always induces topological defects at 
long scales, leading to short range (exponential) correla¬ 
tions in the CDW. Without disorder, the model reduces 
to a pure XY model and it allows a Berezinsky-Kosterlitz- 
Thouless (BKT) phase transition [88, 89] from a low tem¬ 
perature phase with power law correlation to a high tem¬ 
perature phase with exponential correlation. With finite 
RFD, the low temperature “phase” readily has exponen¬ 
tial correlation, and the residue of the phase transition 
is a crossover [81], below which only short range static 
CDW order exists. Now we return to the full model con¬ 
taining PDW2, CDW and RFD, which has been studied 
before in Ref. [65]. From the results of random field XY 
model, we readily conclude that the CDW has exponen¬ 
tial correlation. This also implies the absence of long 
range PDW order, since long range PDW would have in¬ 
duced CDW with long range correlation, contradicting 
the fact that the system with RFD only has short range 
CDW. Besides, we also expect a crossover temperature 
for the short range static CDW at low temperature. To 
this end, we will review the topological defects [45, 90-92] 
arising in the PDW2 system (section 11.1) and the RFD 
induced crossover (section 11.2). We also present numer¬ 
ical evidence to show the presence of these topological 
defects in the PDW2 system with finite RFD in section 
11.3 and lay the foundation for later numerical studies of 
the four-component PDW systems. 

We first consider PDW2 state and its order parameter, 
which generally is 

A(i?) = A+Q(i?)e+'«-^ + A_Q(i?)e-'«-^ , (1) 

where A±q{R) G C. In particular, if A±g(R) = Ao/2 
are equal and homogeneous in space, then the PDW order 
A(R) = Aq cos Q • R shows modulation of wave vector Q 
in real space. This state induces a CDW order at wave 
vector 2Q, such that P 2 q{R) ^ Alg(R)A+g(R). 

In order to discuss the influence of random field dis¬ 
orders on the PDW and CDW state, we write down the 
Ginzburg-Landau (GL) functional (density) for PDW2 
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[46, 65]: 

Fi^ = +a^ |A|^ + |A|" + 73^ |A+q|^ [A.qI^ 

+«i|a,A|2 + «(^|a,A|2 , (2) 

where A = (A+g, A_g)^. In general, the system could 
admit nematicity and would differ. The other part 
of the GL functional involves the CDW and is given by 

FP = -H*P2Q + C.c. + < \d,p2Q\^ + \dyP2Q\^ , (3) 

where H{R) = hrmsh{R)e'^^^^"> is the local complex ran¬ 
dom field, h{R) are Gaussian-distributed random vari¬ 
ables with mean 0 and standard deviation 1, r]{R) G 
[0,27r) are uniformly distributed random phases, and 
hrms ^ 0 controls the (root mean square) random field 
disorder strength. The random field term models ran¬ 
domly distributed non-magnetic charged impurities that 
pin down the GDW modulations. Notice that the RFD 
couples to the GDW, which refiects the charge modula¬ 
tions of the PDW, but not directly to the PDW. In the 
following, we assume that the PDW always induces a 
CDW p{R) ^ P2q cos 2 Q • R such that 

P2Q = AlgA+g . (4) 

In this case, the total GL functional can be expressed in 
terms of A±g only as 

F2 = F^ + FP (5) 

F'’ = -F*A*_qA+q+c.c. 

+< |4 (A*_qA+q) 1^ + kI \dy (A*_qA+q) f . 

We note that 73^ controls whether A±g components com¬ 
pete with each other. In this section, we shall consider 
7^ < 0 such that two A±g components coexist and in¬ 
duce a 2(5-CDW, even in the absence of random field 
disorder H{r) = 0. The case 7^ > 0 will be discussed in 
section III. 


II. 1. Topological defects 

Here we briefiy discuss the possible topological defects 
[41, 45, 90, 91] in the PDW2 system (in the absence of 
RFD). Assuming 7^ < 0, then we have |A+g| ^ |A_g|. 
We thus write A±g(R) = where we assume 

the amplitude fiuctuations are small compared to those 
of the phases 0±q. The PDW order is now 

A(R) = Aoe*'^ cos {Q ’ R^ p) , (6) 

where the superconducting phase is'd = | (^+Q + ^-o) 
and the PDW smectic phase p = \ (^+Q ~ ^-q)- Deep 
in the PDW2 phase in 2D, the amplitude fiuctuations are 
negligible, the system can be effectively described by an 
anisotropic XY model in terms of the P and p [45, 90-92]. 
The PDW2 system possesses three types of topological 


defects labeled by {uy^rid)'. (1) pure SG 27r-vortices with 
(±1, 0), (2) 27r-dislocations with (0, ±1), and 7r-vortex-7r- 
dislocation defect, or simply half-vortex, with (±^,±|). 
The labels (uy^nd) are related to the topological singu¬ 
larities (topological “charges”) in'd(R) and p{R), defined 
by 

d£' V'd = 27Tny^ 

di'Vp = 27rnd . (7) 

Due to the single-value conditions in the original ^±g(r) 
fields, Uy^d can take half-integer or integer values and 
it leads to three types of defects in the system. More¬ 
over, starting with a PDW2 at low temperature, owing 
to (thermal) proliferation of different defects, the sys¬ 
tem restores broken symmetries and yields various phases 
[45, 91, 92]. A pure GDW phase is accessed when the SC 
27r-vortices are proliferated and a charge-4e SC phase is 
obtained through the proliferation of 27r-dislocations. To 
attain a nematic phase one needs to proliferate the half¬ 
vortices. We note that since the PDW order modulates 
half as many as its induced CDW, the PDW half-vortex, 
manifested as a 27r-dislocation in the CDW, is accompa¬ 
nied by a half SC fiux 4>o/2, which can be experimentally 
verified [65]. 

II.2. PDW2 with RFD 

It is well known that the random field disorder de¬ 
stroys CDW long range order for dimension d < 4 
[62, 63]. Consider the RFD term —i7*AlgA+g + c.c. oc 
— cos{p — 2 p) in Eq. (3), the CDW smectic phase 2 p 
tends to align with the RFD phase 7. Treating the ran¬ 
dom field as a perturbation, the energy gain due to dis¬ 
order potential ^ dominates over the elastic energy 
cost of adjusting to disorder ^ As a result, long 

range CDW is inhibited beyond a length scale known 
as the Larkin length. This argument implies the absence 
of long range PDW order if the CDW is induced by the 
PDW. 

This picture only considers continuous elastic defor¬ 
mation of the uniform state, but not topological defects, 
which are non-perturbative in nature. Indeed, via nu¬ 
merical studies [81, 82], it is shown that these topologi¬ 
cal defects always exist for any RFD strength in 2D. In 
the absence of these defects, the CDW correlation de¬ 
cays in power law, admitting a genuine phase transition 
to a high temperature phase with exponentially decaying 
correlation. In 2D, however, at a length scale 
the equilibrium system is always unstable to a prolifera¬ 
tion of the static topological defects and thus the CDW 
correlation decays exponentially for any RFD strength. 

In 2D, these CDW topological defects manifest as 27r- 
dislocations. Re-tracking back to the PDW, it means 
a PDW half-vortex exists exactly at the CDW 27r- 
dislocation [65]. Proliferation of these half-vortices de- 
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FIG. 3. (Color online) Monte Carlo results for PDW 2 sys¬ 
tem in F 2 , lattice [Eq. ( 8 )]. The phases of the average orders 
(a) (A+q), (b) (A_q), and (c) (p 2 q) are shown. Here we 
choose 73 < 0 such that we can consider the phases and their 
fluctuations only and assume constant order magnitudes. In 
(a) and (b), several defect pairs with opposite windings are 
present in (A±q). Due to the fact that the CDW is induced 
by the PDW, these defect pairs also manifest in the {p 2 q) 
[see (c)]. As an example, panel (d) shows a 2ti CDW smec¬ 
tic phase winding around a counterclockwise path of a CDW 
27r-dislocation (see the insert) located near R = (53,60). Pa¬ 
rameters (in meV): = 5, ^ = 1, /irms = 4, < 0 

and 7 i^ > 0. System size: 96 x 96 with periodic boundary 
condition. Simulated annealing is performed from a high tem¬ 
perature to the final temperature at 10 K. 


stroys the long range PDW and lead to a nematic phase. 
These topological defects are pinned and static (over 
thermal average) at low temperature, while they can dif¬ 
fuse driven by thermal fluctuations at higher tempera¬ 
ture. The presence of these defects results in short range 
correlations in PDW and CDW [65]. Instead of a genuine 
phase transition, we thus expect a crossover temperature 
Tco^ below which short range static charge order can 
be observed, and above which most of the charge order 
is destroyed by the defect diffusion (except domains that 
are pinned strongly by the RFD). Here only an argument 
is given, we shall provide numerical evidence in section 
III. 


II.3. Numerics: Monte Carlo study of PDW2 with 

RFD 


In order to gain more intuition on the interplay of 
PDW, CDW and RFD in a PDW2 system, we perform 
Monte Carlo study on the corresponding lattice GL func¬ 
tional [65], 


E2, lattice — ^"2 ^lattice + ^ 2 , lattice 


( 8 ) 


with 


_ 

-^ 2 ,lattice 

+ 5 ^ ^ |A,Q(i )|2 + ^ |A,Q(i)|^ 

is i \ s 

+7^^l^+Q(i)l^|A-Q(i)|^ 

- 2 K^^Re [AsQ(i + 7 A*Q(i)] 

is 

-2Kf ^ Re [AsQ(i + y)A*Q{i)] 

is 

= 

2 ,lattice 
2 /Arms ^ ^ 

^ Re [A_Q(i + x)A\Q{i + x)AlQ(i)A+Q(i)] 

i 

'F + y)^-Q(i)A+Q(i)] 

+2 (k'i' + <) ^ |A_Q(i)A+Q(i)|^ 

on a square lattice. Here + 2 k^ and 

the summation is over site i and sign s = ±. x and y 
are shifts by one lattice constant in the corresponding 
directions in real space. The random fields are 

Gaussian-distributed at each site. If we choose 73 ^ < 0, 
we can safely set the amplitudes |A±q| as constants for 
the system at low enough temperature and only concern 
about the phase fields 0±q{i). This essentially reduces 
to a random field XY model [65], on which we perform 
the Monte Carlo (MC) simulations with Metropolis al¬ 
gorithm of single-site phase rotations. The initial PDW 
phases 0 ±q are initialized randomly at a high enough 
temperature, which is then slowly lowered to the desired 
one. This simulated annealing process lets the system 
to be better equilibrated at the final temperature. At 
the final temperature, we measure the thermal averaged 
phases of A±q and p 2 Q = A^gA+g at each site. For 
simplicity, we choose = k±. 

In Fig. 3, the results are shown for a typical RFD con¬ 
figuration. Due to their single-value conditions, we show 
the 0±Q fields instead of the and cp fields. The choice 
of parameters deserves some comments. Restricting to 
7 ^ < 0 (i.e. coexisting A±g), the phase diagram for 
^" 2 ,lattice turns out to be simple. In Ref. [65] it is found 
that there are only two phases: the fully disordered phase 
with the presence of all types of topological defects, and 
the PDW glass phase, in which no CDW dislocation of 
27r-winding is allowed (only 47 r-defects are permitted). 
The parameters in Fig. 3 correspond to the fully disor¬ 
dered phase due to the presence of 27r-defects. The PDW 
glass phase is interesting in its own right, but we shall 
only consider the fully disordered phase here. Back to 


h{i)e AlQ(i)A+Q(i) 
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Fig. 3, several observations are in order. First, the CDW 
phase 2(^(i) forms “domains.” This is consistent with 
our previous discussion that a few sites bare strong ran¬ 
dom fields and the CDW phases near these sites would 
align along with these strong random fields under the 
influence of the stiffness terms . Next is the exis¬ 
tence of 27r-dislocations (27r-windings) in the CDW phase 
[Fig. 3(c) and (d)]. We note that these defects occur near 
the domain boundary since doing so the system can lower 
its energy better than pure continuous elastic deforma¬ 
tion [81]. These dislocations come in pairs of opposite 
windings, which is expected in the system with periodic 
boundary condition and zero total winding. 

We note that this simulated method can be utilized to 
find the ground state of the system, but what we have 
found might not be the true ground state and it is not our 
main object to do so. Indeed, the system shows glassy 
behavior that it can be trapped in some local energy min¬ 
imum. This can be understood in terms of the phase re¬ 
lation 2 (p = 0-^Q — 0-Q that although the CDW disloca¬ 
tions are largely determined by the RFD configuration, a 
CDW 27r-dislocation can manifest itself as a 27r-vortex in 
either PDW component A±q [see Fig. 3(a) and (b)]. This 
“selection” depends on the initial configurations and the 
simulated annealing process (for instance, the tempera¬ 
ture intervals during cooling and the MC steps for each 
temperature). Different selection would give a slightly 
different total energy, but it is understood that in 2D, 
the interplay of PDW, CDW and RFD always produces 
a short range PDW and CDW state with the presence of 
topological defects. 


III. PHASE-SEPARATED PAIR NEMATIC 
(PSPN) STATE 

In this section, motivated by the STM and REXS ex¬ 
perimental results, we introduce a four-component PDW 
(PDW4) state with phase separation in directional com¬ 
ponents, which is dubbed as a phase-separated pair 
nematic (PSPN) state. Experimentally, STM results 
[19, 56] reveal that the presence of short range CDW 
within the pseudogap regime. These CDW orders are di¬ 
rectional and have wave vectors 2Qx and 2Qy. It is also 
shown that these two CDW components compete with 
each other in real space and form a domain structure, 
or phase separation pattern. This picture is also sup¬ 
ported by the REXS [13] showing that the results are 
more consistent with stripy CDW domains. Viewing this 
in terms of PDW (section II), within each domain, the di¬ 
rectional CDW with 2Qa {a = x^y) is actually induced by 
a PDW2 with components and this subsystem can 

be described by the PDW2 model Eq. ( 8 ). However, for 
the global system with A^, and Ay domains (A^ denotes 
A±g^ components), we have to generalize to a PDW4 
model [40, 41, 43-48]. Moreover, in order to reproduce 
phase separation pattern, we shall consider strong com¬ 
petition in the directional components A^ [ 7 ^ > 0 , see 




FIG. 4 . (Color online) (a) Schematic density plot of F4,local 
[Eq. (12)] as a function of |A±q| for 7 ^ > hrms > 0. (b) Plot 
of F4,local as a function of © defined in (a) along the white 
dashed line with |Asq|^ = ~2^ ~ F4,local has two 

minima separated by an energy barrier si. At either mini¬ 
mum, the induced CDW p 2 Q and LC order I are both finite. 


Eq. (10) below]. 

On the other hand, TRSB order is detected in the 
pseudogap by PKR [25, 54] and PND [27-29]. The ex¬ 
periments seem to support that the loop current (LC) 
order [49] is responsible for the TRSB order. A natural 
question would be whether we can explain this TRSB 
order within the PDW4 model. Interestingly, Agterberg 
et al. [48] proposed that the LC order can arise from a 
PDW theory as a secondary composite order. However, 
it is readily discovered that within the PDW4 model the 
LC order is incompatible with the induced CDW in a 
PDW state with spatially homogeneous components Ag. 
One way to resolve this issue is to introduce an eight- 
component PDW model, as done by Agterberg etal. [48], 
while another feasible way, as will be demonstrated be¬ 
low, is to include RED and to relax the assumption that 
the state has homogeneous PDW components. 

In the PDW4 model, another relevant parameter is 7 ^ 
[see Eq. (10) below] controlling whether the A+g^ and 
A_g^ components compete. If the components A±q^ are 
assumed to be spatially independent, from the known re¬ 
sults on the possible PDW4 ground states [41, 48], we are 
constrained to choose 7 ^ < 0 such that A+g^ and A_g^ 
can coexist. Otherwise, no CDW p 2 Q^ = Alg^A+g^ 
can be induced, contradicting the experimental results 
discussed above. By introducing disorders in the form of 
RED, we break the translational invariance and allow the 
PDW components A±q^ to be spatially inhomogeneous. 
We will demonstrate that in this setting, the parameter 
range 7 ^ > 0 (for competing A+g^ and A_g^) should 
also be taken into account. Inclusion of RED and the 
parameter choice of 7 ^, 7 ^ > 0 lead to two important 
consequences: (i) generating a phase separation in P 2 g^,^ 
(induced by A^^y) that is consistent with that observed 
in the experiments and (ii) resolving the incompatibility 
issue of CDW and LC orders within the PDW4 model. 
Eurthermore, we show that in the PDW4 model there is 
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also a crossover temperature Too < 7^* to account for the 
short range static CDW within the pseudogap. For tem¬ 
perature Tco < T < T"", the induced CDW is dynami¬ 
cally fluctuating. The simulation result also suggests that 
the PDW superconductivity is absent owing to strong SC 
phase fluctuations, complying with non-superconducting 
pseudogap. These results constitute the properties of the 
PSPN state. 

We consider the PDW4 system on a square lattice. 
We shall assume that the two new components A±g^ are 
related to the PDW2 components A±g with Q = by a 
7 r/ 2 -rotation. By construction, we require that Qx ^ Qy 
The order parameter is given by 


A(i) = ^A,Q„(i)e-‘3“-. (9) 

as 

We write down the lattice PDW4 GL functional [46] by 
generalizing the PDW2 case. 


7^4,lattice — ^4,lattice + ^4,lattice 


( 10 ) 


with 

rA _ 

■^4,lattice 

- 2 Kf ^Re [A,q„ (i + ail) A.y^ (i)] 

ias 

-2k^ Re [Asq„ (i + a_L) A*q^ (i)] 

ias 


+a^5]|A,Q„ (i)|2 

ias 

i \ as / 

ia 


+ 27 ^ ^ Re [a+q^ (i)A_Q^ (i)A+Q^ (i) 


PP 

4,lattice 


- 2 /irmsy)Re ha{i)e 

ia 

-2kI Y. Re [A_q„ ■ A;q^ (i + ail)A*_ q^ ■ A+q„ (i)] 

ia 

-2< ^ Re [A_q„ • A;q^ (i + a^) A*_q^ • A+q, (i); 

ia 

+2 ^ |A_g^(i)A+g^(i)|^ , 


where we have implicitly replaced p2Qa Alg^A+g^ 
and Alg^- A+g^(i) denotes Alg^(i)A+g^(i) for brevity. 
Here the summation is over a = and + 


The notation ay = a and a± = y (x) if 
a = X (y). The rest is similar as in the PDW2 case, 
except that now we have four PDW components at :^Qx,y 
and allow coupling ( 7 ^ term) between the directional 
components. Notice that each 2Qa-CDW is coupled to a 
different random field disorder 

In the absence of RFD hrms = 0, the ground state of 
the PDW4 functional is mainly controlled by 72 and 73 . 
As in the PDW 2 case, the sign of 7 ^ controls whether two 
components A+g^ and A_g^ (a = x, y) coexist ( 7 ^ < 0 ) 
or compete ( 73 ^ > 0) with each other. While the sign 
of 7 ^ dictates the coexistence ( 7 ^ < 0 ) or competi¬ 
tion ( 7 ^ > 0 ) of \Axf and \Ayf components, where 
|Aa| = l^sQal • literature, the main concern 

is on the latter. (There is also a similar debate [93, 94] on 
pure CDW model) In Ref. [43-45], the authors propose a 
PDW theory with 7 ^ > 0 and 7 ^ < 0 such that in a clean 
2D system, the ground state is a (unidirectional) striped 
PDW with either A^, or Ay component. On the other 
hand. Ref. [47] suggests a bidirectional PDW state with 
7 ^, 73 ^ < 0 that the Ax^y components coexist. Experi¬ 
mentally, X-ray scattering measurements [13, 15] reveal 
the CDW with only wave-vectors 2Qx and 2Qy but not 
A2Qx ± 2Qy. If the PDW, as well as the induced CDW, 
is bidirectional, then it is likely that the A2Qx ± 2Qy 
CDW components would be observed [93] (see, however. 
Ref. [47]). Here we propose another parameter choice 
72 ^, 73 ^ > 0 such that all four components A^g^ com¬ 
pete with each other and constitute the PSPN state. We 
shall show that it can generate, besides the CDW orders, 
a LC order [49] accounting for the polar Kerr effect and 
the intra-unit cell time-reversal symmetry breaking order 
observed respectively in the PKR and PND experiments. 

The (phenomenological) LC order emerged from PDW 
order A±g^ is defined as [48] 

la = \A+Qf-\A.Qf. ( 11 ) 

This order parameter is translational invariant (hence 
giving rise to an lUC order observed in PND), odd under 
either time reversal symmetry and parity, and invariant 
under their product. (We note that from a pure CDW 
order, one can also construct a similar composite TRSB 
order [95, 96]) Since this LC order is originated from the 
PDW, so its characteristic temperature is the same as 
T*, which is indeed observed in the PND [27-29]. At 
first sight, in the PDW4 model this order la is incompat¬ 
ible with the CDW order P 2 Qa' If we restrict ourselves 
to spatially independent order, then a finite la would re¬ 
quire 7 ^ < 0 , which leads to one of the A+g^ and A_g^ 
components to vanish and eventually the induced CDW 
order P 2 Qa — 0 - This is the reason why the authors in 
Ref. [48] turn to an eight-component PDW theory. For 
details of the eight-component theory, the readers shall 
refer to the reference. Here we will demonstrate that if 
we relax the assumption of spatial homogeneity of the or¬ 
der parameters and introduce RFD, finite orders in both 
la and P 2 q^ are actually compatible within the PDW4 
model. 
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Before diving into the numerics, we shall first consider 
a simplified model to understand the interplay of (short 
range) CDW and LC order with the RFD. The simplified 
“on-site” model is given by 

^4,local = (E (12) 

+7^ I^+qI^ - 2ft.rms |A_qA+(5| 

such that we can ignore the phase degrees of freedom. 
Here we have already assumed 7 ^ > 0 and only one pair 
of A±q exists locally. We further assume that the CDW 
smectic phase is already aligned with the local pinning 
field such that the RFD term reduces to the present form. 
Now we notice that although 7 ^ > 0 would suppress one 
of the A±q components, but the RFD term {hrms > 0) 
induces a finite CDW order ~ |A_qA+( 5 | such that both 
7 ^ 0- This competition, provided with the suitable 
parameters (see Fig. 4), can result in a local state with 
both non-vanishing p 2 Q and L We shall briefly remark 
that the effect of RFD on the eight-component PDW the¬ 
ory [48], in which the incompatibility issue of CDW and 
LC order is solved by introducing extra competing PDW 
components. If RFD is introduced, then using a similar 
argument as above, it is likely that extra CDW compo¬ 
nents other than those of wave vectors 2Qx and 2Qy can 
be induced. While currently there is no experimental 
evidence showing these extra CDW components. On the 
other hand, in Fig. 4 we notice an interesting feature that 
the two minima are separated by an energy barrier si. 
One may define a temperature scale Tk ^ (^/(i))i as the 
spatial average of the energy barriers, which characterizes 
the stiffness of the state to (thermally activated) relative 
amplitude fluctuations between |A+g^| and |A_g^|. For 
T < Tk, the energy barrier hinders the relative ampli¬ 
tude fluctuations and this results in finite LC order. For 
T > Tk^ the thermal fluctuations overcome the energy 
barriers and the system now consists of fluctuating LC or¬ 
der domains with different signs leading to a thermal av¬ 
erage la) ^ ^ for TRSB order. This argument might 
explain the disparity between the temperature scales for 
TRSB order measured in PKR and PND, in which the 
one measured in PKR is strictly lower than that in PND 
(see Introduction). The elaborated study however will 
be pursued elsewhere. 

Now back to the entire system described by F^^iatUce 
[Eq. 10] with RFD , we focus on those sites 

with particularly strong pinning strengths. At one of 
these sites, to gain the pinning energy, the CDW smec¬ 
tic phase would align with the RFD phase 7 (see section 
11.2). The PDW and CDW stiffness terms would induce 
a finite region of approximately the same smectic phase 
around the strong disorder site. Now if we are restricted 
to the parameter choice discussed above and take the 
directional PDW into account, we would then expect a 
phase separated domain pattern, in which the directional 
PDW, together with the induced directional CDW and 
LC orders, presents due to the “nucleation” around these 


strong disorders. The phase separated orders discussed 
above can be regarded as imposing phase separation pat¬ 
terns on the PDW2 orders. This indicates that the cor¬ 
relations of the PDW and CDW orders are bounded by 
those in PDW2, and thus the resultant state also has 
PDW and CDW with the short range exponential cor¬ 
relations. Analog to the PDW 2 case, we also expect a 
crossover temperature for the short range static CDW in 
the state. Owing to the phase separation pattern and 
the local nematicity originated from PDW superconduc¬ 
tivity, we shall call this the phase-separated pair nematic 
(PSPN) state. Below we will demonstrate that the PSPN 
state can be resulted from the PDW4 GL functional with 
RFD [Eq. 10] via MC simulation study. 


III.l. Numerics: Monte Carlo study of PSPN 

To better illustrate our points, we perform the MC 
simulations similar to the PDW2 case. We consider the 
case of 72^3 > 0 as we discussed above and choose 
and such that the average jA^g^l ^ 1 . We divide 
the MC study into two stages. In the first stage, in or¬ 
der to capture the effects of 72^3 on the resultant state, 
we allow updates (Metropolis algorithm, single site) on 
both the amplitudes and phases of the PDW order pa¬ 
rameters. We start with a random phase configuration 
with |Asg^| = 1 and perform the simulated annealing 
process to let the system to better equilibrate at the fi¬ 
nal temperature. For later comparison to experimental 
data, the simulation is performed with a typical RFD 
configuration. We then obtain a state at low tempera¬ 
ture, resembling the ground state for a given set of pa¬ 
rameters. In the second stage, we are interested in the 
effect of thermal fluctuations on the orders, in particular 
the phase fluctuations induced by proliferation of topo¬ 
logical defects (see section II.l). Analog to the studies of 
BKT transition [ 88 , 89] (see also random field XY model 
[65, 81]), we shall assume that the system is deep in the 
state obtained in stage one. We thus utilize the aver¬ 
aged amplitudes as input (i.e. the amplitudes are kept 
fixed during the MC) and only perform MC simulations 
with phase updates (Metropolis algorithm, single site) for 
various temperatures. To be consistent, and used 
in stage one should be chosen to suppress the the total 
PDW amplitude fluctuations but still permit sufficient 
acceptance rates for sampling states with various rela¬ 
tive PDW amplitudes. Besides, directional CDW with 
wave vectors 2 Qx,y requires that 7 ^ is dominated over 
7 ^, while in order to induce finite orders in both CDW 
and LC order, 7 ^ and hrms should be comparable (see 
Fig. 4). We shall set 7 ^ = 0 such that the phases of the 
directional PDW are not coupled [97]. 

In Fig. 5, we show the simulation results in stage one 
for a typical RFD configuration for PDW4. Several ob¬ 
servations are followed. Owing to a large 7 ^ > 0, the two 
components Ax^y compete with each other and it results 
in phase separation for Ax^y^ as well as (Recall 
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FIG. 5. (Color online) Monte Carlo results of the PSPN state in the presence of RFD for the averaged order parameters in (a) x- 
and (b) ^-directions, using CL functional F 4 , lattice Eq. (10). We plot the spatial dependent (I and II) PDW components (A±q^), 
(III) CDW {p2Qa) and (IV) LC order {la). The PDW and CDW orders are phase separated in the directional components. 
Also the PDW and CDW phases arg(AsQ^) and arg(p2Qa) form a domain pattern, similar to the PDW2 case. Moreover, finite 
CDW and LC orders are induced in the system. In order to capture the effects of 7^)3 > 0, we allow updates on both amplitudes 
and phases of the PDW order parameters A^q^. Parameters (in meV): — —100, 7 ^^ = 25, 7 ^ = 20, 7 ^ = 5, 7 ^ = 0, 

= 10, /^ll ^ = 1, /irms = 2 . Systcm size: 96 x 96 with periodic boundary condition. Simulated annealing is performed from 
a high temperature to the final temperature at 10 K. 


that Aa denotes directional PDW components A±q^) 
Within each domain of A^^y and p 2 Q^y^ the phases 
(argA±Q^ and argp2Qa) have similar domain structures 
as in the PDW2 case, except that there are regions of van¬ 
ishing amplitude in which the phases are not well-defined. 
There are still topological defects in the domains, but the 
number is not well constrained due to the regions of van¬ 
ishing amplitude. Experimentally, CDW dislocations of 
27r windings are indeed observed in STM [19] and they 
could trap the half SC fluxes for a PDW state (see sec¬ 
tion II. 1), provided that the SC phase is strongly pinned. 
One can probe these half fluxes by SQUID microscopy 
and test the proposed PDW theory. Next, we notice both 
p2Qa la have finite magnitude locally, in accord with 


our discussion on the simplified model Eq. (12). Also, 
as discussed, the PDW and CDW orders are expected 
to be short ranged and decay exponentially [65, 81] (see 
Fig. 1 in Supplementary Material [98] (SM) for numeri¬ 
cal evidence and more discussions), due to the random 
phases and presence of topological defects. The LC or¬ 
der la are inhomogeneous and have random signs. This 
selection of the local sign is related to the initial config¬ 
uration and simulated annealing process. We also note 
that the amplitudes lA^g^l are inhomogeneous, to ac¬ 
commodate both non-vanishing p2Qa ^a- Again, the 
simulation might be trapped at some local energy mini¬ 
mum and unable to attain the true ground state, but the 
main purpose is to illustrate the effect of RFD and the 
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FIG. 6. (Color online) Melting of short range PDW and 
CDW in the PSPN state in the presence of RFD by per¬ 
forming MC simulations on GL functional F 4 ,lattice [Eq. (10)]. 
As a function of temperature, the averaged local PDW and 
CDW orders, defined as ^ l^sQa(i)P and 

-local _ ^ \p2Qa (01 5 arc plottcd, where N is total number 

of sites. Notice that the plotted values are rescaled by their 
maximal possible values. In order to study the effect of the 
phase fluctuations, the amplitudes of the PDW components 
|Asq^| are taken from the averaged results in Fig. 5 as in¬ 
put and kept constant during the simulations, while only the 
phases arg Asq^ are updated. We set the threshold po = 0.2 
[97] for the local CDW order to define the crossover tempera¬ 
ture Tco for the short range static CDW orders. Parameters 
are the same as in Fig. 5. 



0 48 96 0 48 96 0 48 96 0 48 96 


FIG. 7. (Color online) Spatial variations of CDW in the PSPN 
state at two temperatures: (a) T = 50K < Tco and (b) 
T = 190 K > Tco- In the same MC simulations of Fig. 6 using 
GL functional Fk,lattice [Eq. (10)], we show the results for the 
averaged directional CDW {p2Qa)- The noticeable difference 
for these two temperatures is the vanishing small magnitudes 
of the averaged CDW at T = 190 K due to the strong thermal 
fluctuations in the smectic phases for T > Tco- 


simultaneous presence of la and p 2 Q^ for our choice of 
7# > 0- 

We then move on to stage two to study the effect of 
the thermal phase fluctuations on the PDW and CDW 
orders. As discussed, owing to the pinning of the RFD, 
the PDW and CDW orders are short ranged, and we ex¬ 
pect a crossover temperature Tco for the static CDW. To 
determine Tco^ we plot the temperature dependent local 
CDW order as well as local order for PDW, 

in Fig. 6. The local orders decay in a continuous fashion 
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FIG. 8. (Color online) Spatial variations of PDW components 
in the PSPN state at two temperatures: (a) T = 50K < Tco 
and (b) T = 190K > Too- In the same MC simulations in 
Fig. 6 using GL functional F 4 ,lattice [Eq. (10)], we measured 
the PDW order parameters (A±q^ ^). At T = 50K, some 
domains of directional PDW are wiped out owing to strong 
fluctuations in SC phases. At T = 190 K, the fluctuations are 
stronger and most of the PDW orders are destroyed. These 
results suggest the absence of PDW superconductivity irre¬ 
spective of Tco- 


as temperature increases. At low temperature, the local 
charge orders with significant amplitude are still present, 
as shown in Fig. 7(a), due to the strong pinning effect of 
the RFD and correlation effect of the stiffness terms . 

At high enough temperature, the focal charge orders are 
almost absent [Fig. 7(b)], as the thermal phase fluctua¬ 
tions have overcome the pinning effect from RFD and 
stiffness terms. We can view this as the thermally acti¬ 
vated diffusion of the topological defects (dislocations) of 
the CDW. Moreover, although the amplitudes are kept 
fixed as the averaged values taken from Fig. 5, we no¬ 
tice that as a manifestation of the strong phase fluctu¬ 
ations, the domains of well-defined CDW smectic phase 
at high temperature [Fig. 7(b)] shrink considerably. We 
shall roughly identify the crossover temperature for the 
short range static CDW Tco = 150 K by setting the 
threshold po = 0-2 [99]. In this sense. Fig. 7(a) and (b) 
actually correspond to the spatial variations of the direc¬ 
tional CDW for T < Tco and T > Tco respectively. 

Next we discuss whether the PDW superconductivity 
can survive the thermal phase fluctuations. Fig. 8 shows 
the MC results for the magnitudes and phases of the aver¬ 
aged PDW orders (Asq^) at T = 50 K and 190 K. At T = 
50 K < Tco [Fig. 8(a)], there are PDW domains with al¬ 
most vanishing amplitudes, in contrast to the same area 
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in Fig. 7(a), in which charge orders of significant ampli¬ 
tude still remain. This difference is because of the strong 
pinning on the CDW smectic phase by the RFD, while 
PDW superconducting phases+ ^-Qa) 
not (here OsQ^ = argA^g^ are the phase fields of the 
PDW components, see also section II. 1). Since the PDW 
smectic phases ip a = — 0-q^) are also pinned in¬ 

directly by the RFD through the induced CDW below 
Tco^ fhis suggests that the fluctuations in the SC phases 
'da in the domains are strong. We note that the phase 
fluctuations are stronger than the PDW2 case due to the 
phase separation pattern [97]. The thermal fluctuations 
can “invade” the orders at the boundaries and destroy 
them more severe than the system without phase separa¬ 
tion. At T = 190 K > Tco [Fig-8(b)], the thermal fluc¬ 
tuations are severe and the PDW local order is mostly 
wiped out, as evident from the vanishing averaged PDW 
order parameters and the considerable shrinkage in the 
PDW domains of well-defined phases. It is well known 
that dynamic vortices in a superconductor dissipate su¬ 
percurrent, leading to a finite resistance (known as the 
flux flow resistance^ see e.g. Ref. [100]). In this regard, 
the results in Fig. 8 hint that the PDW superconduc¬ 
tivity is absent owing to strong thermal fluctuations in 
the superconducting phase irrespective of the crossover 
temperature Tco and appear to be consistent with the 
known fact that the pseudogap is non-superconducting. 

We conclude this section by comparing these results 
with experiments. First, RIXS [8, 9, 11] and STM [17] 
reveal traces of the crossover temperature Tco for static 
charge order inside the pseudogap. RIXS uses hard X- 
rays to measure the static order of a given structure. Sim¬ 
ilarly STM measures the static electronic orders. T^o’s 
found in these probes show a similar trend that is peaked 
at some doping. Comparing with our results, we identify 
the order to be the short range CDW. REXS measure¬ 
ments however use soft X-rays, that is sensitive to fluctu¬ 
ating orders, and it can detect the fluctuating CDW up 
to a higher temperature T*, at which the PDW orders 
vanish. Next, a smectic modulation measure 0^(r, e) 
in STM [18, 19, 56] signifies relative strength of the lo¬ 
cal smectic (directional) electronic orders. As discussed 
at the section overview, these results point to a phase- 
separated directional PDW state. We found that using 
our model of PDW4 and RFD together with suitable pa¬ 
rameters, we can reproduce the phase separation pattern 
for the directional CDW (see Fig. 3 in SM [98] for more 
details). We note that REXS result [13] also supports 
such a phase separation picture for the directional CDW. 


IV. CONSTRAINT FROM CDW WITH 
PREDOMINANTLY D-WAVE FORM FACTOR 

Before the discussion of the spectral functions of the 
PDW system, we shall first figure out the pairing sym¬ 
metry of the PDW. Both STM [20, 56] and REXS [15] 
experiments strongly indicate that the CDW order within 


the pseudogap has a predominantly d-wave form factor. 
If we insist that the CDW order is only originated from 
the corresponding PDW order, then we need to under¬ 
stand how this CDW form factor would constrain the 
pairing symmetry of the PDW. 

Here we propose that the PDW has a 5 '± id-symmetry. 
Eirst we write down the definition of a generic SC order 
parameter at sites and Vj 

A(ri,r2) = ^ 

Q \ k / 

(13) 

in which we allow multiple components with different 
modulation wave vectors Q, which determine the mod¬ 
ulation in the average coordinate R = | (ri -hr2), and 
their (orbital) pairing symmetries involving relative co¬ 
ordinate r = ri — r 2 are described by the form factor 
AQ{k). We can similarly define the order parameter for 
charge modulations [20] by a simple replacement A ^ p. 
We note that the definition above only describe the order 
parameters of the PDW components Ag spatially uni¬ 
form in the average coordinate (i.e. Ag is independent of 
R). If we want to include (slow) spatial variation in Ag, 
we can add a spatial varying factor by the substitution 
Ag(r) ^ Ag(R)/g(r) with the form factor encoded in 
/g(r). In this section, we shall consider spatially uniform 
case only. 

The pairing symmetry factor for TQx is assumed to be 

A+Q^(fc) = A,(fc2 + kl) + - kl) , (14) 

where As^d G M , 0 G [0, 27r). It is clear that the first (sec¬ 
ond) term is 5-wave (d-wave). Here Ag^d and (j) will be 
constrained by the CDW dominating d-wave form fac¬ 
tor shown in the experiments. Time reversal symme¬ 
try requires that A_g(/c) = A;^g(—/c). By the relation 
P2Qa{'^) = ^-g^(^)^+Qa^be CDW form factor can 
be shown as (see Appendix A) 

P2QAk) = {Al + e^^’^Al/2)ikl + klf 

^ ^ ^ 

s-wave 

+i2d^A,A,){kl+kl){kl-kl) 

' -V-' 

d-wave 

+{e^^^Al/2){kt-6klkl + kl) (15) 

'-V-' 

g-wave 

up to an overall factor. The STM experiment [20, 56] 
reveals that the 5-wave component is much smaller than 
the d-wave one, while the REXS experiment [15] shows 
that a pure d-wave CDW cannot fully reproduce the ex¬ 
perimental data. Inspired by these experimental facts, 
we simply choose the 5-wave component A^ -h \e^'^^A^ 
to vanish and leaving only the d-wave and ^-wave com¬ 
ponents. This particular choice requires = Ti and 
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Ad/As = V 2 . Explicitly we take 


^sQai^) = ^0 


{kl + kl) + isV2{kl 



(16) 


where Aq G C. Here and A_g^ transform to 

each other under time reversal symmetry and we simply 
choose A^g^ GG AsQy related by a 7r/2-rotation. In the 
experiments, the d-wave form factor also indicates that 
the charge orders are located at the oxygen site. There¬ 
fore, in the (single band) model on a square lattice, we 
shall consider the (nearest neighbor) bond order between 
the lattice sites (with Cu atoms) to mimic this experi¬ 
mental fact, and we restrict ourselves to the s' Aid-wave 
pairing symmetry 

^sQai^) _ (^cos /Ca,+COS /c^)+25^/2(cos /Ca,—COS /c^) . (17) 
Ao 

A detailed derivation of the induced CDW form factor for 
Eq. (17) is given in the Appendix A, in which we show 
that the above consideration remains valid in the lattice 
scenario. 


V. QUASIPARTICLE SPECTRAL FUNCTIONS 
AND ARPES 

In this section, we study the spectral functions of PDW 
systems and compare the results with the ARPES exper¬ 
iments. Many exotic properties are observed in ARPES. 
Early studies have revealed that the Eermi surface in the 
pseudogap forms the so-called Eermi arcs [12, 51], discon¬ 
nected segments of gapless quasiparticle excitations near 
the nodes. Near the antinodes, the spectrum is gapped, 
while the quasiparticle peak is very broad or even ill- 
defined [2]. Recently, it is found [54, 55] that by com¬ 
paring the data above and below T*, the quasiparticle 
spectra have the feature of kp-kc misalignment near the 
antinodes, where kp is the Eermi momentum for T > T"" 
and kc is the back-bending momentum for T < . This 

feature cannot be explained by conventional SC pairing, 
otherwise we would expect kp = kc- Moreover, away 
from the antinode, the spectral gap is closed from below 
[55, 101, 102]. While from a calculation [47] for a model 
with pure CDW, an order observed in STM and REXS, 
shows that the spectral gap is closed from above and it 
fails to capture this ARPES feature. On the other hand, 
PDW models [47, 48] can successfully capture the kp- 
kc misalignment and the feature of the gap closing from 
below, suggesting that PDW plays an important part in 
the pseudogap phenomenology. 

To compute the spectral functions, we introduce a lat¬ 
tice hamiltonian with PDW pairings [Eq. (18)] taking 
into account the spatial variations in the PDW order 
parameters, pairing symmetry and thermal (phase) fluc¬ 
tuations. In order to study separately the effects of the 
RED induced random phases, phase separation, and the 
thermal fluctuations, we will consider four PDW systems: 
(i) uniform PDW2 with s' ± id-wave pairing symmetry. 


(ii) s' ± id-wave PDW2 with random phases induced by 
RED, (iii) s' ± id-wave PDW4 with phase separation in 
pairing amplitudes Ax^y but no random phase, and lastly 
(iv) PSPN state with thermal fluctuations in the PDW 
phases Osq^ , originated from a 5 ' ± id-wave PDW4 state. 
Eor the former three systems, it is sufficient to compute 
the spectral functions by feeding a static PDW configu¬ 
ration into the hamiltonian [Eq. (21)]. But for last one, 
it is needed to incorporate disorder average and thermal 
fluctuations through a sequence of MC generated PDW 
configurations [Eq. (22)]. 

Actually there are similar studies on the effect of ther¬ 
mal fluctuations on the spectral features. In Ref. [103], 
the authors studied the thermal averaged spectral func¬ 
tions for cuprates within a spin liquid scenario, where 
the order parameter fluctuates between a d-wave super¬ 
conductor and a nearby staggered-flux state related by 
SU(2) rotations. While in Ref. [104], the authors studied 
the fluctuating d-wave superconductor scenario, in which 
the SC phases are strongly disordered due to vortex pro¬ 
liferation. These studies also appear to reproduce the 
Eermi arc feature observed in ARPES. 

We first write down the real space hamiltonian H = 
Hq + with 

Hq — ~ + h.C. — fi <^icr^icr (18) 

iSa icT 

HA='^A*{i + a, i) (ci+afCu - Ci+a^Cit) + h.C., 

ia 

where A (i -h a, i) is the nearest neighbor singlet (bond) 
pairing with a = x^y. The parameters are taken from 
Ref. [12] that (with lattice constant uq = 1) the nearest 
neighbor hopping ti = 0.4 eV, the next nearest neigh¬ 
bor hopping t 2 /ti = —0.2 and the third nearest neigh¬ 
bor hopping ts/ti = 0.05. The PDW order parame¬ 
ters in momentum space A±q^ (k) describe pairing with 
total momentum AQa and the pairing centers are at 
AQa/^. In order to explain the ARPES results, we fol¬ 
low Ref. [47] and choose the pairing centers at the Bril- 
louin zone (BZ) boundary that ±Qx/‘2 = ±(7r,(3o/2) 
and AQy/2 = ^ {Qq/2^p) related by a 7r/2-rotation (see 
Eig. 1). In this study, we take the chemical potential 
/i/ti ~ —0.758 (corresponding hole doping ^ 0.117) such 
that the electronic band intersects with the BZ bound¬ 
ary exactly at the aforementioned pairing centers with 
Qo = tt/A. The resultant PDW has a period of 8ao and 
corresponding CDW with period 4ao. As (hole) doping 
increases, the corresponding CDW modulation 2Qq de¬ 
creases by choice [47] and this trend is in accord with the 
results observed in X-ray experiments [10, 12, 14]. 

In order to incorporate spatial non-uniform pairing or¬ 
der parameters, we generalize Eq. (13) to 

A(ri, r 2 ) = A y] AQ{R)fQ{r)e ^^'^, (19) 

Q 

where R = ^(ri + r 2 ) is the average coordinate, r = 
ri — r 2 is the relative coordinate, /g(r) describes the 
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pairing symmetry in real space, given by [Eq. (17)] 


fsQair) 


1 + isV2, r = ±x 
1 — is^/2 , r = ±y 


( 20 ) 


Here the spatial dependent pairing order parameters 
Aq{R) G C will be taken from the results of the MC 
simulations with spatial mean magnitude {\Aq{R)\)ji = 
1. We shall assume that the orders are sufficiently 
smooth and we can take the interpolation value Ag (R) = 
^ (Ag(ri) + Ag(r 2 )) as the bond order parameters [105]. 
Finally, an overall scaling A is added to control the “gap” 
size. In the following, we compute the quasiparticle spec¬ 
tral functions in momentum space corresponding to the 
hamiltonian Eq. (18) according to 


^{A} (^: ^) = ——Im (^k 


1 


uo ■ 


i0+ - H. 


{A} 


k) , (21) 


provided that the pairing order parameters {A} = 
{A^g^} are given. Here \k) = 

single-particle wave vector with momentum k in the first 
BZ. 


V.l. PDW2 

We start with the uniform s' ± id-wave PDW2 system 
[Eq. (8), 7 ^ < 0]. In particular, we choose to add the 
A±Qy PDW orders and the results of selected momen¬ 
tum line scans [see Fig. 1] are shown in Fig. 9(a). More 
line scans are shown in Fig. 5 of SM [98]. We first notice 
that the Fermi surface is gapped at both the antinodal 
regions {kx,y = tt scans), despite that we have only added 
the PDW pairings centered at AQyj2. The spectra near 
the nodal regions {kx,y = f scans) are gapless and re¬ 
semble the bare spectra without the PDW pairing, ex¬ 
cept the existence of some gap structures away from the 
Fermi energy. These results constitute the Fermi arcs 
observed in the ARPES experiments [12, 51]. Within 
the PDW picture, these Fermi surface segments are due 
to poor pairing condition near the Fermi surface. How¬ 
ever, the “gap,” defined as the energy gap at the back- 
bending momentum kc below the chemical potential, is 
not monotonically decreasing from the antinodes [com¬ 
pare Fig. 5(a-IV) and (a-V) in the SM [98]]. This is how¬ 
ever inconsistent with the ARPES [52, 53]. Next, the 
band structures near two antinodes are different. Near 
AQx/2 [Fig. 9(a-HI)], the spectrum resembles closely a 
conventional band structure gapped by an order param¬ 
eter such that the lower band edge at ^ — 0.2eV is left 
below the Fermi energy. Near AQyl2 [Fig. 9(a-I)], where 
the PDW order is centered, we notice the so-called kp- 
kc misalignment of the band gap, where kp is the Fermi 
momentum of gapless bare bands (white dash lines) at 
high enough temperature. This misalignment observed 
in the ARPES experiments imposes a strong constraint 
on theories. Although this uniform PDW2 state can ex¬ 
plain several important features observed in ARPES, but 


more satisfactory results could be obtained for the sys¬ 
tems we will consider below. 

Next we turn to the s' ± id-wave PDW2 system with 
RFD induced random phases. We set the amplitude 
\A±Qy{R)\ to be constant and use the averaged phases 
in Fig. 3 obtained by MC simulation as input to compute 
the spectral functions. The results with the same line 
scans are shown in Fig. 9(b). The spectra are similar to 
those in the uniform PDW2 state, except a noticeable 
difference that the coherent peaks are more broadened 
near the pairing centers ±(5^/2 than at the other antin¬ 
odal region AQxl2 [compare Fig. 9(b-I) and (b-HI)]. Near 
the nodal regions, the spectra are largely unaltered by 
the random phases, except some minor spreading of the 
spectral weights. So the main effect of random phases is 
broadening of the coherent peaks near the pairing cen¬ 
ters and the band structure induced by the PDW orders 
is largely preserved. 


V.2. phase-separated PDW4 

We now discuss about the s' ± id-wave PDW4 state 
with the Aa,^^-phase separation [Eq. (10), > 0]. 

Firstly, we notice that the PDW2 spectra in Fig. 9, es¬ 
pecially near the antinodes, do not preserve the fourfold 
rotational symmetry, while this symmetry appears to be 
unspoiled in the ARPES results. In the following, we 
show that the RFD induced phase separation pattern for 
PDW4 can help to resolve this discrepancy, namely ap¬ 
proximately preserving the fourfold rotational symmetry 
while still maintaining all the spectral features obtained 
in PDW2. 

Here the PDW orders exists at all four momenta 
^Qx,y‘ We take the averaged Ax^y amplitudes with phase 
separation in Fig. 5 as input to calculate the spectral 
functions. Also focusing on the effect of phase sepa¬ 
ration, we set the PDW phases as spatial independent 
(i.e. constant). The results of line scans kx,y = tt are 
shown in Fig. 10(a) and more scans are presented in Fig. 6 
of the SM [98]. PDW orders now exist at all four pair¬ 
ing centers AQx^yj^^ the phase separated orders result in 
spectrum broadening near both antinodal regions. The 
effect of amplitude fluctuations in PDW4 is similar to 
that of phase fluctuations in PDW2 above in terms of 
spectrum broadening, except it affects both antinodal 
regions. The nodal Fermi surface is again gapless with 
minor peak broadening [Fig. 6(a-VH) and (a-XIV) in SM 
[98]]. Next we observe that the spectra appear as su¬ 
perposition of those from a PDW2 with A±g^ orders 
and another one with A±g^ orders (see Fig. 7 in SM [98] 
for detailed discussion) and thus approximately preserve 
the fourfold rotational symmetry. This is an intriguing 
result since the RFD induces directional PDW domains 
breaking the local fourfold rotational symmetry. We note 
that the spectra at kx,y = tt [Fig. 10(a)] actually agree 
with the ARPES results [54, 55] better than the PDW2’s 
[Fig. 9(a-I)], in which the lower band edge near ^ —0.2 eV 
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(a) PDW2 without rand, phase 


(b) PDW2 with rand, phase 



FIG. 9. (Color online) Spectral functions at two line scans (see Fig. 1) in BZ. Panel (a) shows the PDW2 system without adding 
random phases and panel (b) shows the PDW2 with static random phases taken from Fig. 3. The PDW orders A±q^ are added 
at two pairing centers The white dash lines represent the bare bands. Both systems show the kp-kc misalignment in 

the ky = TV plots as seen in subplots (a-I) and (b-I), where kp is the Fermi momentum of the bare band and kc is the back 
bending momentum. But the spectra in the system with random phases are more broadened. We take Qo = Tr/d, A = 0.125 eV, 
ti = 0.4 eV, t 2 /ti = —0.2, ts/ti = 0.05, and ji/ti ^ —0.758. System size: 96 x 96 with periodic boundary condition. More line 
scans are provided in the Supplementary Material [98]. 


(a) phase separated PDW4 
without thermal fluct. 


(b) PSPN 

T = 5()K 


(c) PSPN 

T = 190K 



(a-l) ky = TV 


(a-ll) kj: = TV 


-0.4 0.0 0.1 -0.4 0.0 0.4 -0.4 


FIG. 10. (Color online) Spectral functions for (a) phase-separated PDW4 without thermal fluctuations, (b) PSPN state with 
thermal fluctuations at T = 50K < Tco, and (c) PSPN state at T = 190 K > Tco are plotted at kx,y = tv. The PDW orders 
/S.±Q^ y are added at four pairing centers In panel (a), we take the averaged PDW amplitude results |A±q^ from 

Fig. 5 as input and the PDW phases are homogeneous and constant. The spectra closely resemble superposition of two PDW2 
spectra with either A±q^ or A±q^ and the fourfold rotational symmetry is approximately preserved. In panels (b) and (c), 
we have performed the thermal and disorder averages at the stated temperatures. We note that the spectrum superposition 
feature revealed in (a) remains. In all cases, the kp-kc misalignment feature is intact, though the broadening is more severe in 
the two high temperature cases. The parameters are the same as in Fig. 9. More line scans are provided in the Supplementary 
Material [98]. 
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is missing. We also note that the dispersion approaches 
Fermi level from below when moving from antinodes to¬ 
wards the nodal regions [see Fig. 6 (a) and Fig. 8 in SM 
[98]], consistent with ARPES results [55, 101 , 102 ]. Near 
the end of the gapless Fermi arc (see Fig. 8 in SM [98]), 
the dispersion is not particle-hole symmetric but bends 
downward and loses spectral weight, qualitatively con¬ 
sistent with the ARPES results by Yang etal. [101, 102 ]. 
The same data, as well as some antinodal features, has 
also been nicely captured by the phenomenological YRZ 
model [ 102 , 106-108]. We have also performed a similar 
calculation on the PDW4 with coexisting ( 7 ^ < 0), 
the spectrum does not have the superposition feature and 
fails to reproduce the ARPES results due to the strong 
interference effect of the Qx and Qy PDW orders. 


V.3. PSPN state with (thermal) phase fluctuations 


We are finally ready to discuss the spectral function 
results for the PSPN state with phase fluctuations. This 
state stems from the s' ± id-wave PDW4 state. In order 
to capture the phase fluctuations of the PDW orders, it 
is insufficient to compute only the static averaged con¬ 
figurations as in above. Actually we need to compute 
the averaged qausiparticle spectral function over thermal 
fluctuations, given by [104] 


{A{k,ui)) 


^{A} 


A{A}{k,Uj)e ^ 4 ,lattice 


( 22 ) 


where for brevity {A} denotes {Agg^}. This averaged 
spectral function {A{k,uj)) is calculated by taking the 
thermal average with {Agg^} configurations generated 
by MC simulations of the GL functional Ei^iattice in 
Eq. (10). We then perform the disorder average over 16 
RED configurations. We carry out the numerical calcu¬ 
lations at two temperatures T = 50K and T = 190 K, 
which are at opposite sides of the charge order crossover 
temperature of Too ^ 150 K. Again, for each RED con¬ 
figuration, we take the amplitudes from the MC results 
with both amplitude and phase updates and simulated 
annealing process as in Eig. 5 as input and then carry 
out the MC simulations on the PDW phases at the de¬ 
sired temperatures. The results for 50 K and 190 K are 
respectively presented in Eig. 10(b) and (c). The spectra 
are smoothened compared with those without phase fluc¬ 
tuations, but the band structures are largely unaltered 
compared with the PDW4 results without phase fluctu¬ 
ations [Eig. 10(a)]. Many important features that we ob¬ 
tained previously remain intact even with the presence of 
the phase fluctuations. Eirstly, the “Eermi arcs” remain. 
The spectra near the antinodes are severely broadened. 
Compared with the PDW4 results without phase fluc¬ 
tuations [Eig. 10(a)], the broadening mainly comes from 
the Aa,^^-phase separation instead of the phase fluctua¬ 
tions. The “gap” decreases monotonically away from the 


antinodes owing to the phase fluctuation induced gap fill¬ 
ings and finally becomes gapless near the nodal regions. 
Again away from antinodes, the spectrum approaches the 
Eermi level from below, in accord with ARPES experi¬ 
ments [55, 101 , 102 ]. Secondly, the band structure su¬ 
perposition and the kp-kc misalignment also survive the 
phase fluctuations. Moreover, if we compare the results 
for 50 K and 190 K, we see that the crossover temper¬ 
ature Tco does not induce a qualitative change on the 
spectrum. The two major features mentioned above still 
survive the depinning of the smectic phases across Tco 
and can be regarded as definitive characteristics of the 
PDW orders for T < T"". These results agree well with 
the ARPES measurements [ 12 , 52-55, 109, 110 ], in which 
also show these features in the pseudogap region below 
T*. (We note that these ARPES features can also be 
captured to some extent by alternative theories [78, 111]) 
Comparison of the result and the ARPES data is given 
in Eig. 9 of the SM [98]. 

By comparing the spectral results for different PDW 
states, we notice two major effects of RED. Eirstly, in 
terms of the phase separation pattern and the random 
phases, it broadens the spectrum, as one would have 
expected from a disordered system. The broadening is 
more severe near the antinodes, which is also observed 
in ARPES data [ 2 ]. But more importantly, RED induces 
the phase separation pattern for the directional orders, 
rendering the spectra approximately obeying the fourfold 
symmetry, while still preserving many experimentally ob¬ 
served features appeared in the directional PDW2 state. 

Lastly we argue that the superposition feature actu¬ 
ally support the proposal that the cuprate SC has stripe 
PDW and CDW orders ( 7 ^^ > 0). Since without PDW 
order, the fermionic spectrum in the CDW phase cannot 
reproduce the ARPES results [47, 55, 101 , 102 ] that the 
spectrum approaches from Eermi level from below, we 
shall confine ourselves to various scenarios of PDW theo¬ 
ries. Eirst, the pure PDW 2 could not reproduce the lower 
band edge of the ARPES spectrum as we discussed above. 
Next, we have also checked that the checkerboard PDW4 
state with 72^3 < 0 could not reproduce the ARPES re¬ 
sults. A third possibility is proposed in Ref. [47] that the 
ARPES results actually constitute from a checkerboard 
PDW4 orders decaying in the momentum space. It is pro¬ 
posed that the PDW symmetry factor fsQ^ (fc) ^ e-^ /« 
decays away from the pairing centers, and it can suc¬ 
cessfully explain the kp-kc misalignment and the Eermi 
arcs. But the spectrum is largely coming from the PDW2 
states at the pairing centers due to the exponential de¬ 
caying SC symmetry factor, so it has missed the ARPES 
feature at the lower band edge shown in Eig. lO(a-I). Nev¬ 
ertheless, the stripe PDW does not necessarily need to 
exist within the same layer. Actually, alternating layers 
of X- and ^-directional stripe PDW [40] would produce a 
similar superposition spectra. But the phase separation 
structure similar to those in Eig. 5 is indeed observed in 
STM [19, 56] and also supported by REXS [13], so it is 
perhaps favored over the alternating layer scenario. It is 
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hard to exhaust all the possibilities, but it is intriguing to 
see the convergence of the STM, REXS (indicating stripe 
CDW, induced by the stripe PDW) and the ARPES (in¬ 
dicating stripe PDW) experimental results towards a uni¬ 
fying theoretical framework. 


VI. CONCLUSIONS 

In conclusion, the phase separated pair nematic 
(PSPN) state derived from a four-component pair density 
wave (PDW) model under the influence of random field 
disorders (RED) captures coherently a number of experi¬ 
mental results in the pseudogap. We consider, to the best 
of our knowledge, a new parameter regime 7^3 > 0 for 
the phenomenology of cuprates, in which all the PDW 
components are competing with each other. The 

new choice of parameter impedes the generalization of 
existing PDW results to the present case, however, we 
show that in terms of a phenomenological GL approach, 
the experimental observed features obtained in the exist¬ 
ing PDW theory are retained if disorders are included. It 
is shown that CDW and a time-reversal symmetry break¬ 
ing (TRSB) order in the form of loop current (LC) are 
induced as secondary composite orders. The inclusion of 
disorders permits the state to be spatially inhomogeneous 
(in terms of the PDW components) and this readily re¬ 
solves an issue concerning the incompatibility of CDW 
and loop current order within a four-component PDW 
model [48]. Nonetheless, the CDW is short ranged under 
the influence of random field disorders, consistent with 
the STM [19, 20, 56] and REXS [13] results, while the LC 
order can account for the TRSB order observed in polar 
Kerr rotation [25, 26] and polarized neutron diffraction 
[27-30] experiments. Eurthermore, by a MC simulation, 
we show that random field disorders can also induce a 
phase separation pattern for the CDW, similar to that 
observed in STM. This PSPN state also appears to cap¬ 
ture other experimental features. It explains a distinct 
temperature scale Tco < for the static (short range) 
CDW (see section III). More importantly, it accounts for 
the same pseudogap temperature T* observed in linear 
resistivity, ARPES (for the anomalous spectral features), 
REXS (for the dynamically fluctuating directional CDW) 
and polarized neutron diffraction (for the lUC TRSB or¬ 
der) experiments. We also argue that the thermal su¬ 
perconducting phase fluctuations lead to finite flux flow 
resistance, resulting in a non-superconducting state re¬ 
gardless of Tco- To test the proposed PDW theory, one 
may probe the half SC flux trapped at the 27r-CDW dis¬ 
location (see section II.I and III.I), provided that the 
local SC phase is strongly pinned. In section IV, we con¬ 
strain the PDW, from the observed CDW with a pre¬ 
dominantly d-wave form factor, to have s' ± id pairing 
symmetries. In section V, we show that a number of 
anomalous features in ARPES, namely Eermi arcs, kp- 
kc misalignment, antinodal gap closing from below are 
retained in the PSPN state, even with thermal fluctua- 


Acronym 

Full phrase 

ARPES 

Angular resolved photoemission 
spectroscopy 

BZ 

Brillouin zone 

CDW 

Charge density wave 

dSC 

d-wave SC 

lUC 

Intra-unit cell 

LC 

Loop current 

MC 

Monte Carlo 

NMR 

Nuclear magnetic resonance 

PDW 

Pair density wave 

PDW2 

Two-component PDW 

PDW4 

Four-component PDW 

PKR 

Polar Kerr rotation 

PND 

Polarized neutron diffraction 

PSPN 

Phase separated pair nematic 

REXS 

Resonant inelastic X-ray scattering 

RIXS 

Resonant elastic X-ray scattering 

RED 

Random field disorder 

SC 

Superconductor / superconductivity / 
superconducting 

SM 

Supplementary material (Ref. [98]) 

STM 

Scanning tunneling microscopy 

TRSB 

Time reversal symmetry breaking 


TABLE I. Acronym list for the phrases commonly used in the 
paper. 


tions. Moreover, the PSPN state is shown to have severe 
coherent peak broadening near the antinodes (but not 
the nodes), which has been a puzzling feature observed 
in ARPES [ 2 ]. The random field disorder induced phase 
separation pattern helps to explain why the ARPES spec¬ 
trum still obey the fourfold rotational symmetry approx¬ 
imately, while locally it has been broken as revealed in 
STM. Erom above, we see that the random field disor¬ 
ders assist in understanding the pseudogap for 7 ^ < 7 ^* 
in terms of the PSPN state by resolving several issues in 
previous studies. 
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Appendix A: CDW form factor induced by s' zb id 
PDW on a square lattice 

From the definitions of PDW Eq. (13) with symmetry 
factor Eq. (14) 

f+Qik) = {kl + kl) + a{kl-kl) 
f.Q{k)=flQ{k) 

(here a G C) and its induced CDW 

P {r, +2Q) = A* (r, -Q) A (r, +Q) , 

we have 

p{r,2Q) oc ^ /+Q(fci)/+Q(^ 2 )e-*'=^''e+*''^’' 

ki k‘2 

P 

and 

p27r 

P2q{p)= deqf+Q{q+p/2)f+Q{q-p/2) 

Jo 

where 9q = arg g, and we define ki = g+1 and /c 2 = g — |. 
After some straightforward calculations, we have 

P2Q{p) = J^[0+a^/2)ipl+plf 

+2a{pl+pl){pl-pl) 

+{ay2){pt-Qplpl+pl)]+0{p^) 

We note that the first term is clearly 5-wave. While the 
second (third) term is d-wave (g-wave) due to the idem 

titles cos26'p = and cos46»p = ■ A less 

rigorous derivation by the product of the basis functions 

P2q(p) ~ (1 + Q!Cos26>p)^ 

= (1 + a^/2) + 2acos2dp + {c? jlP) cos46*p 


would also give the same result. 

Similarly, in the lattice case, we have symmetry factor 
Eq. (17) 

/(^) = (cos kx + cos ky) + i\/2 (cos kx — cos ky) . 
Then 


p{p, 2Q) = 

q 

+ 1 — 2y/2i^ cos kiy cos k 2 y 

+3 cos kix cos k 2 y + 3 cos kiy cos k 2 x • 

Then we substitute /ci=g+|,/c 2 =g—| and evaluate 
the integrals ^ dqx dqy, we have 

p(p, 2 ( 5 ) = +A|f^ {cospx - cospy) 

{cospx + cosp^) , 

where and A^^^ = —27r^. In the continu¬ 

ous case, after setting the s'-wave component to zero, we 
have 

p2Qik) = 2V2iAl{ki - kl) - A^ik^ - 6klkl + k^). 

The ratio of the d-wave to ^-wave component is —2\/2i, 
matches that of the lattice case A^^^/A^g^ = —2^/2i. We 
note that although the induced CDW’s form factor has a 
“^-wave” component, it manifests as 5 '-wave ^ cospx + 
cosp^ due to the choice of nearest neighbor bond order 
and the square lattice. 


(-1 + 2 V 2 ij COS kix cos k2x 
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I. SHORT RANGE DIRECTIONAL CDW AND CHARGE-4E SUPERCONDUCTING CORRELATIONS 

IN PSPN STATE 

In this section, we discuss the correlations of directional CDW and charge-4e superconducting order in the phase- 
separated pair nematic (PSPN) state. In section I.l, we shall briefly review these correlations and the corresponding 
phases in the two-component PDW system. Next in section 1.2, we generalize the correlations to the four-component 
PDW system and provide numerical evidence for the short range correlations as argued in the main text. In section 
1.3, we discuss the rationale behind the parameter choice 74 ^ = 0 and its implications to the PSPN state. 


I.l. Brief review on two-component PDW model 


Here we briefly summarize the possible states derived from a two-component PDW (PDW 2 ) system with components 
We define the composite orders, namely the CDW 


P2Q{r) = AlQ(r)A+Q(r) 


and the charge-4e superconductivity 


A 4 e(r) = A_( 3 (r)A+Q(r) 


( 1 ) 

( 2 ) 


As briefly discussed in the section II. 1 in the main text, there are fou r possible state s induced by proliferating different 
topological defects provided that the local PDW order |A(r)| = is non-vanishing (see for example 

Ref. [? ? ] for more details). The proliferation of SC 27r-vortices destroys the superconductivity, resulting in a (pure) 
CDW state. In this CDW state, the directional CDW correlations 


Cp,2Q{r - r') = (p2Q(r)p2Q(r')) (3) 

are long ranged and the charge-4e SC correlation 

CA, 4 e(r - r') = (A2g(r)A4e(r')) (4) 

is short ranged (exponential decay). We note that the short range charge-4e SC correlation is induced by the prolif¬ 
eration of SC 27r-vortices, as in the Berezinsky- Kosterlitz-Thouless (BKT) transition. We summarize the states, the 
corresponding topological defects proliferated, and the correlations in Table I as the ease of reference. In particular, 
we shall focus on the nematic state with short range Cp^ 2 q{t) and CA, 4 e(^) is induced by proliferation of the half¬ 
vortices, which would manifest as CDW 27r-dislocations in p 2 Q and charge-4e SC 27r-vortices in A 4 e (see Table II for 
more details). 


1.2. Four-component PDW model and PSPN state 

In the four-component PDW system, we generalize the PDW 2 case to define the composite orders, namely the 
directional CDW 


/32Q„(r) = AlQjr)A+Q„(r) 


( 5 ) 



2 


{a = x^y) and charge-4e superconductivity 

A4e{r)= • 

a=x,y a 


( 6 ) 


We note that since only phase separated directional PDW components = y^|A+q^P + |A_g^p are considered 
here, so we simply add the directional charge-4e order parameters 

^4e W = ^-Qa (^)^ + Qa (^) . (7) 

to give A 4 e as above. Similar to the PDW 2 case, there are four possible states induced by proliferating different 
topological defects provided that the local PDW order |A(r)| = non-vanishing. 

For the PDW4 system, we similarly define the directional CDW correlations 

Cp, 2 Q^ {r - r') = {{p* 2 Q^ ir)p 2 Q^ (r'))) (8) 

and the charge-4e SC correlation 

C'A,4e(r-r') = ((A24r)A4e(r'))). (9) 

Here ((•)) denotes thermal and disorder averages. Depending on the behaviors of these correlations we can also have 
four possible states as in the PDW2 case. In particular, we shall focus on the nematic state with short range Cp^2Qa{'^) 
and CA, 4 e(^) is induced by proliferation of the half-vortices, which would manifest as CDW 27r-dislocations in p2Qa 
and charge-4e SC 27r-vortices in A^g. 

Next, we describe the MC process to compute the correlations under the influence of the PDW phase fluctuations. 
For each RFD configuration, we first perform a simulated annealing down to T = 10 K with both amplitude and phase 
updates (as in stage one described in section III.l in the main text) to obtain the averaged PDW orders. We then use 
the averaged PDW orders to perform a MC simulations and measure various correlations shown in Fig. 1 with only 
PDW phase updates (as in stage two described in section III .2 in the main text). Finally, we average the results over 
a number of RFD configurations to obtain the disorder average. 

The results on the Cp^2Qa{'^) C'A, 4 e(^) a^re shown in Fig. 1. It is clear that both correlations decay exponentially 

and are short ranged, resulting in a nematic state (at least locally within each directional PDW domain). These 
correlations are actually induced by the half-vortices in the directional PDW, which manifest as CDW 27r-dislocations 
in p 2 Qa [see Fig. 5(a-III) and (b-III) in the main text] and charge-4e SC 27r-vortices in Ate (see Fig. 2). 


1.3. Rationale for choosing 7r = 0 

We finally discuss the implication of the fact that the pseudogap is non-superconducting on the parameter 7 ^ in 
the GL functional Eq. (10) of the main text. For the ease of reference, we repeat the 74 ^ term here 


2^4 (^)A+q^ (^)A-q^ (r) 


( 10 ) 


We note that if we identify the directional charge-4e SC components A 4 g(r) = A_Q^(r)A+g^(r), then the 74 ^ term 
can be regarded as the Josephson junction coupling between the local Afg and A^g orders. We shall consider the case 
7 ^ is finite and negative. Then we expect that the directional charge-4e SC orders A^g would attain phase coherence 
at low enough temperature, even if their amplitudes are phase separated as shown in Fig. 2 . Notice that the 

total charge-4e SC order Eq. ( 6 ) is finite everywhere, since we require the finiteness in the PDW orders. As in the case 
of BKT transition, the phase coherence in A 4 e then implies that the system is superconducting below T < ^ 7 ^ 

(i.e. a charge-4e superconductor, see e.g. Ref. [? ]). Experimentally, the pseudogap region in the temperature range 
Tc < T < T* is not superconducting, where Tc is the transition temperature for the d-wave SC. This fact suggests 
that 7 ^ is small, such that the charge-4e superconductor induced by the phase separated PDW orders discussed above 
is hidden inside the d-wave SC dome, i.e. T 4 e < Tc. In the main text, we simply set to zero In this simplified case, 
the phase coherence between directional charge-4e SC A^g cannot build up owing to phase separation pattern (see 
Eig. 2 ) and we expect the PSPN state to be non-superconducting. The above discussion provides a rationale for such 
a parameter choice. 
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Type of topological defect proliferated 

Gharge-4e SG correlation 

GDW correlation 

PDW 

none 

Long range 

Long range 

GDW 

(1,0): SG 27r-vortex 

Short range 

Long range 

Gharge-4e SG 

(0,1): 27r-dislocation 

Long range 

Short range 

Nematic 

(|,|): half vortex 

Short range 

Short range 


TABLE L List of four states accessible in a two-component PDW (PDW2) system by proliferation of different types of topological 
defects labeled by (jiv^nd) in Eq. (7) in the main text. 


PDW 

(|,|): half vortex 

(1,0): SG 27r-vortex 

(0,1): 27r-dislocation 

GDW 

GDW 27r-dislocation 

- 

GDW 47r-dislocation 

Gharge-4e SG 

Gharge-4e SG 27r-vortex 

Gharge-4e SG 47r-vortex 

- 


TABLE IL List of topological defects in a two-component PDW (PDW2) system labeled by (riv^nd) in Eq. (7) in the main text. 
In the second and the third rows, we list the topological defects in terms of the induced CDW and charge-4e superconductivity 
manifested by the corresponding defects in the PDW. 
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FIG. 1. Spatial correlations of (a,b) directional CDW p2Qa — ^-Qa^+Qa (c) charge-4e superconductivity 

A 4 e = ^a=x y ^-Qa^+Qa f^c PSPN statc at T = 10K. The results are averaged over 80 RFD configurations (disorder 
average). The blue dots are the mean values of the correlations and the error bars are for one standard deviation. The red 
dashed lines are the exponential fits for each set of data. The results show that the PSPN state has both short range CDW and 
charge-4e SC correlations and it is nematic (locally). See the discussions in section I above for the meaning of the correlations, 
corresponding states and the MC process. 



FIG. 2. Spatial variations of averaged directional charge-4e superconducting orders (A^g) measured at T = 10 K, where 
A 4 e(r) = A_q^( r)A+Q^(r). It is obvious that phase separation also occurs in the directional components A^^^, inherited from 
the directional PDW orders. It is also evident that the phases of directional A^g are not coupled due to the choice of = 0. 
From the averaged orders at low temperature, we see that static charge-4e SC 27r-vortices is present owing to the pinning 
of the PDW indirectly through the coupling between RFD and CDW. This result can actually be inferred from the averaged 
PDW orders (Asq^) shown in Fig. 5 of the main text. These static vortices lead to the short range charge-4e SC correlation 
CA, 4 e{r) shown in Fig. 1. The MC process is the same as those in Fig. 1, except that only one RFD configuration is computed 
here. 
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FIG. 3. Comparison of the Monte Carlo (MC) result for the CDW (a) \{p2Qx)\ and (b) |(p 2 Q^)| [taken from Fig. 5(a-III) 
and (b-III) in the main text] and (c) the smectic modulation measure 0?(e) = Z {e)~^ [ReZ{Sy,e) — ReZ{Sx,e)] in the STM 
experiments [1, 2], where we shall identify Sx{y) — ‘^Qx{y)- Here Of{e) signified the relative strength of the local smectic 
(directional) CDW order in the x- and y-directions. We note that the MC simulation largely reproduces the STM results 
that show phase separation of 2Qx- and 2Q2/-CDW’s. Note: In Ref. [1] and [2], the authors defined another quantity, the 
intra-unit-cell electronic nematicity: 0^(e) = Z{e)~^ [ReZ((5^,e) — ReZ((5^,e)], where Q'^ — (27r/ao,0) and Qy — (0,27r/ao). 
This quantity does not measure any smectic modulation since Q'x^y is just the reciprocal lattice vectors (not to be confused 
with our Qx,y for the PDW wave vectors). 



FIG. 4. Line scans in blue for the quasiparticle spectral functions in Fig. 5 and 6. Those in red for Fig. 8. 
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(a) PDW2 without rand, phase 


(b) PDW2 with rand, phase 



FIG. 5. The figure is an extended version for Fig. 9 in the main text and more line scans (see Fig. 4) are presented. Panel (a) 
shows the PDW2 system without random phases and panel (b) shows the PDW2 with static random phases taken from Fig. 3 
in the main text. 
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(a) PDW4 without thermal fluct. 

(a-l) ky = TT (a-VIII) = TT 


(b) PSPN at T = 5()K 

(b-l) ky = n (b-VIII) = TT 


(c) PSPN at T = 190K 

(C-l) ky = TT (C-VIII) = TT 
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FIG. 6. The figure is an extended version for Fig. 10 in the main text and more line scans (see Fig. 4) are presented. Spectral 
functions for (a) phase-separated PDW4 without thermal fluctuations, (b) PSPN state with thermal fluctuations at T = 50 K < 
Tco, and (c) PSPN state at T = 190 K > Tco- We note the so-called “Fermi arcs,” that the spectra near the anti-nodes are 
gapped. While moving towards the nodal regions the gap diminishes monotonically and the spectrum eventually becomes 
gapless and coincides with the bare band. 
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FIG. 7. Spectral function superposition feature for the phase-separated PDW4 state. Panels (a) and (b) are the sepctral 
functions at the BZ boundaries along ky = n and kx = t: respectively. Panel (a) [(b)] is taken from Fig. 5(b-I) [Fig. 5(b- 
VIII)]. Panels (c) is for phase-separated PDW4 state, taken from Fig. 6(a-I). We notice panel (c) resembles a spectral function 
superposed from (a) and (b), by noticing that the various features presented in them also appear in (c). We have marked these 
features in black arrows for the ease of reference. We thus conclude that though the PDW4 state is phase separated with Qx 
and Qy domains, its spectral functions resemble superposition spectrum from those of Qx- and (3y-PDW2 states. 
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FIG. 8. Diagonal scans of spectral functions for the phase-separated PDW4 state without thermal fluctuations [see Fig. 10(a) 
in the main text or Fig. 6(a) above]. The numbering of the line scans corresponds to the red lines in Fig. 4. Each line scan is 
in the direction away from the F point and the /c = 0 point coincides with the reduced Brillouin zone boundary. The spectra 
should be compared with those in the ARPES by Yang etal. [3, 4]. Eirstly, the gap near the antinode is closed from below, 
a feature that can be reproduced with PDW but not CDW order [5]. Secondly, the dispersion near the the end of the Eermi 
arc is not particle-hole symmetric but bends downward and loses spectral weight. Both features appear to be in qualitative 
agreement with the ARPES results. Inclusion of random phases and thermal fluctuations seem not affecting the dispersion 
except some broadening as we have seen in Eig. 6, so we expect the results remain valid as well for the PSPN state. 
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FIG. 9. Comparison of the spectral functions for PSPN state and the ARPES results. Panel (a) is taken from Fig. 9(b) and 
shows the spectral functions along different cuts in the BZ from antinode (left) to node (right). The white dashed lines indicate 
the bare bands in different BZ cuts. Panel (b) is taken from ARPES experiment [6] (see also Ref. [7]) and shows the peaks of 
the quasiparticle excitations. Blue and green circles are the peaks for T < Tc and red circles (bare band) are for T > T*. The 
ARPES results consist of denser cuts near the antinodal region, so we cannot directly compare each pair of subgraphs. But we 
can compare the trend and the details of the spectral functions. First, in both cases, the spectrum is gapped near the antinode 
and becomes gapless near the node, forming the so-called Fermi arc. Next, both show the presence of the kp-kc misalignment 
near the antinode. Also, the gap near the antinode is closed by the spectral weight moving towards the Fermi level from below 
when moving away from the antinode [5, 8], which is also captured by the PSPN state in (a). 































